/*
**  File: AquaHow.c
**
**  This file is part of the AQUA program suite, which is being
**  developed in the NMR Spectroscopy Department, Bijvoet Center for 
**  Biomolecular Research, Utrecht. 
**
**
**  AUTHOR(S) : Ton Rullmann
**  VERSION   : 1.5.0
**  DATE      : Dec 18, 1997
*/

#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include "Compiler.h"
#include "General.h"

#include "AquaMacros.h"
#include "AquaTypes.h"

/****************  function prototypes  ****************/

#include "AquaFuncts.h"
#include "AquaFuncts_io.h"
#include "AquaFuncts_dr.h"
#include "AquaStrucset.h"
#include "AquaFiles.h"
#include "MenuFuncts.h"
#include "Functs.h"

#include "AquaHow.h"


/****************  globals  ****************/

extern Boolean  iaf, batch;                           /* import from Functs */
extern Card     msg_string;                           /* import from Functs */
extern int      terseness;                            /* import from Functs */
char            *program_name = "AquaHow";            /* export to Functs */
char            *program_vers = aqua_version;         /* export to Functs */

extern Filnam   dinf_filnamG, dinf_filnamNormal;      /* import from AquaFiles */

int             exit_flag = 0;

enum { opt_skip, opt_sort };
enum { sort_off, sort_only, sort_class };

/* BEGIN_EXTH */

/****************  main  ****************/
/*
**  Determine HOW to analyze a strucset.
**
**  Special exit values:
**    10: file has 0 restraints
**    20: not all restraints were handled because some atoms were not found 
**    30: stereo list is empty
**    40: not all stereo definitions could be made (only for 'stereo assign')
*/

int main( int argc, char *argv[] )
{
    char         *fct_name = "AquaHow main";
    Strucset     theSet;
    int          item, item2;
    int          s_count, dist_count, skip_count = 0, a_count, count;
    Qtyp         qtype;
    DistdataList inf_list;
    Distinf      *dist_inf = NULL;
    Distinfl     *dist_infl;
    Distrst      *dist_rst = NULL;
    DistinfTyp   ftype = inf_normal;
    DRList       DR_list;
    char         *DR_type = NULL;
    Card         DR_ident, lib_string;
    Filnam       dinf_filnam, drst_filnam;
    Filnam       filnam_s, filnam_i;
    Boolean      err;
    Boolean      done     = false, done2   = false;
    Boolean      done_mol = false, done_dr = false;
    Boolean      average  = true,  satur   = true;
    Menu         *oldmenu;
    char         temp_value[4];
    
    enum { item_read, item_mol, item_noe, item_ss, item_hb, item_dis, item_make,
           item_map2aq, item_stereo, item_show, item_out, item_sort, item_select };
    static MenuItem menulist[] =
    {
        { item_EXEC,   "exec",    none, on },
        { item_SET,    "set",     none, on },
        { item_SETS,   "sets",    none, on },
        { item_read,   "readinf", none, on },
        { item_mol,    "mol",     none, on },
        { item_noe,    "noe",     none, off }, 
        { item_ss,     "ss",      none, off }, 
        { item_hb,     "hb",      none, off }, 
        { item_dis,    "dis",     none, off }, 
        { item_map2aq, "maptoaq", none, off }, 
        { item_make,   "make",    none, off }, 
        { item_stereo, "stereo",  none, off }, 
        { item_show,   "show",    none, off }, 
        { item_out,    "out",     none, off }, 
        { item_sort,   "sort",    none, off }, 
        { item_select, "select",  none, off }, 
        { item_QUIT,   "quit",    none, on } 
    };
    static MenuItem menulist_stereo[] =
    {
	{ item_assign,  "assign",   none, on }, 
	{ item_average, "average",  none, on }, 
	{ item_sum,     "sum",      none, on }, 
	{ item_satur,   "satur",    none, on }
    };
    static MenuItem menulist_rstinp[] =
    {
        { item_NONE,   "none",      none, on }, 
        { item_aqua,   "aqua",      none, on }, 
        { item_disman, "disman",    none, on }, 
        { item_disgeo, "disgeo",    none, on }, 
        { item_diana,  "diana",     none, on }, 
        { item_xplor,  "xplor",     none, on }, 
        { item_xplori, "xplor-inv", none, on }, 
        { item_biosym, "biosym",    none, on }, 
        { item_whatif, "whatif",    none, on }, 
        { item_local,  "local",     none, on }, 
        { item_user,   "user",      none, on }, 
        { item_copy,   "copy",      none, on }
    };                      /* same as in AquaWhat, AquaPseudo and GetLIBName */
    static MenuItem menulist_sort[] =
    {
        { sort_off,    "none",   none, on }, 
        { sort_only,   "only",   none, on }, 
        { sort_class,  "class",  none, on }
    };
    static MenuItem menulist_rest_select[] =
    {
        { restsel_both,  "both",  none, on }, 
        { restsel_one,   "one",   none, on }
    };
    static MenuItem optionlist[] =
    {
        { opt_batch,   "batch",       logical, on }, 
        { opt_terseness, "terseness",   integer, on }, 
        { opt_prtaq,   "print",       logical, on }, 
        { opt_prtres,  "printres",    logical, on }, 
        { opt_prtlib,  "printlib",    logical, on }, 
        { opt_prtrst,  "printrestr",  logical, on }, 
        { opt_rstinp,  "restrinp",    oneword, on }, 
        { opt_atmlib,  "atomlib",     string,  on }, 
        { opt_sterlib, "stereolib",   string,  on }, 
        { opt_oldmap,  "oldmap",      logical, on }, 
        { opt_shift,   "shift",       integer, on }, 
        { opt_skip,    "skip",        logical, on }, 
        { opt_guess,   "guess",       logical, on }, 
        { opt_multi,   "multiple",    logical, on }, 
        { opt_sort,    "sort",        oneword, on }, 
        { opt_rest_sel,"rest-select", oneword, on }, 
        { opt_find2,   "find2",       logical, on }, 
        { item_END,    "end",         oneword, on }
    };

    static Menu    menu =
    {
        "Main menu", "command", "Give command:", 
        sizeof(menulist)/sizeof(MenuItem)
    };
    static Menu menu_stereo =
    {
        "Stereo menu", "keyword", "Give keyword:", 
        sizeof(menulist_stereo)/sizeof(MenuItem)
    };
    static Menu menu_rstinp =
    {
        "Restrinp menu", "name", "Give restraint input type:", 
        sizeof(menulist_rstinp)/sizeof(MenuItem)
    };
    static Menu menu_sort =
    {
        "Sort menu", "keyword", "Sort only, or classify as well:", 
        sizeof(menulist_sort)/sizeof(MenuItem)
    };
    static Menu menu_rest_select =
    {
        "Rest-select menu", "keyword", "How to select restraints:", 
        sizeof(menulist_rest_select)/sizeof(MenuItem)
    };
    static Menu options =
    {
        "Options menu", "option", "Give option:", 
        sizeof(optionlist)/sizeof(MenuItem)
    }; 
    

/* initialize */

    ProgramInit( argc, argv );
    ProjectInit();
    StrucsetInit( &theSet, "theSet", true );
    SetMenu( &menu, menulist );
    SetOptions( &options, optionlist );
    SetSubMenu( &options, opt_rstinp, &menu_rstinp, menulist_rstinp );
    SetSubMenu( &options, opt_sort, &menu_sort, menulist_sort );
    SetSubMenu( &options, opt_rest_sel, &menu_rest_select, menulist_rest_select );
    ModOption( opt_rstinp, "none" );
    ModOption( opt_skip, "on" );
    ModOption( opt_guess, "on" );
    ModOption( opt_multi, "on" );
    ModOption( opt_sort, "none" );
    sprintf( temp_value, "%d", terseness );
    ModOption( opt_terseness, temp_value );    
    inf_list.count = 0;
    InitDRList( &DR_list ); 

/* the main menu loop */

    while ( !done )
    {
        item = ReadMenuItem();
        switch ( item )
        {
            case item_EXEC:
		DoExec();
                break;
            case item_SET:
		ReadOption();
		batch = OptionTrue( opt_batch );
                if ( !Check_General_Options () ) {
	            printf( "Error in main: general options set are not okay\n" );
	            MError( );
                }
                terseness   = GetOptionIVal ( opt_terseness );                
                break;
            case item_SETS:
                while ( ReadOption() )
                    batch = OptionTrue( opt_batch );
                break;
            case item_mol:
                if ( ReadAquaSet( &theSet, -1 ) ) /* need all models for mapping.....
		                                     just value of struc_count would do */
                {
                    DisableMenuItem( item_mol );
                    DisableMenuItem( item_read );
                    EnableMenuItem( item_noe );
                    EnableMenuItem( item_ss );
                    EnableMenuItem( item_hb );
                    EnableMenuItem( item_dis );
                    EnableMenuItem( item_show );
		    if ( done_dr )
		    {
			EnableMenuItem( item_make );
		    }
		    done_mol = true;
                }
                break;
            case item_noe:
            case item_ss:
            case item_hb:
            case item_dis:
		DR_type = GetMenuItemKey( item );
                if ( ProcessDRInput( &DR_list, DR_type, DR_ident ) )
		    if ( MapDRInput( &DR_list, theSet.molecule, theSet.struc_count ) )
		    {
			DisableMenuItem( item_noe );
			DisableMenuItem( item_ss );
			DisableMenuItem( item_hb );
			DisableMenuItem( item_dis );
			DisableMenuItem( item_read );
			EnableMenuItem( item_map2aq );
			EnableMenuItem( item_out );
			EnableMenuItem( item_select );
			EnableMenuItem( item_sort );
			if ( done_mol )
			{
			    EnableMenuItem( item_make );
			}
			done_dr = true;
		    }
                break;
	    case item_map2aq:
		oldmenu = SetOneMenu( "Give library (or - to end):" );
		done = false;
		while ( !done2 )
		{
		    MReadWord( lib_string );
		    if ( *lib_string == '-' || strlen( lib_string ) == 0 )
			break;
		    if ( !MapAtomsInDRestraints( &DR_list, ext_to_int, 
		              lib_string, theSet.molecule, theSet.struc_count ) )
		    {
			ERROR( "map failed - destroying the restraint list" );
			free( DR_list.list );               /* destroy the list */
			InitDRList( &DR_list ); 
			ResetMenu( oldmenu );
			EnableMenuItem( item_noe );
			EnableMenuItem( item_ss );
			EnableMenuItem( item_hb );
			EnableMenuItem( item_dis );
			if ( !done_mol )
			    EnableMenuItem( item_read );
			DisableMenuItem( item_map2aq );
			DisableMenuItem( item_out );
			DisableMenuItem( item_select );
			DisableMenuItem( item_sort );
		    }
		    if ( OptionTrue( opt_prtrst ) )
			PrintDRList( &DR_list );
		    done2 = true;
		}
		ResetMenu( oldmenu );
		break;
	    case item_read:
		*dinf_filnam = '\0';
		ProcessDistinfFile( &inf_list, &dist_count, &s_count, &qtype, theSet.filnam, drst_filnam,  
			            DR_ident, false, &ftype, dinf_filnam );
		if ( ftype == inf_normal ) {
		    dist_inf = (Distinf *) (inf_list.fst_desc->data);
                } else {
		    FATAL( "Wrong file/data type - Code error" );
                }
		if ( IsEmptyStr( drst_filnam ) )
		    FATAL( "Name of distance-restraints file not found on distance-info file" );
		DR_type = QtypSel( qtype );
		dist_rst = ProcessDistrstFile( &count, &qtype, filnam_s, filnam_i, DR_ident,  
			                       false, drst_filnam );
		if ( NOTEQ( theSet.filnam, filnam_s ) ||
		     NOTEQ( dinf_filnam, filnam_i )   ||
		     dist_count != count )
		    FATAL( "Inconsistent information on restraints and info files" );
		if ( !ProcessStrucsetFile( &theSet, -1 ) )
		    exit( 1 );
		strcpy( dinf_filnamG, dinf_filnam ); /* save the name of the (normal) info file */
		DisableMenuItem( item_mol );
		DisableMenuItem( item_read );
		DisableMenuItem( item_noe );
		DisableMenuItem( item_ss );
		DisableMenuItem( item_hb );
		DisableMenuItem( item_dis );
		DisableMenuItem( item_map2aq );
		EnableMenuItem( item_show );
		EnableMenuItem( item_stereo );
		MakeDRList( &DR_list, dist_inf, dist_rst, dist_count, theSet.molecule->atoms );
		skip_count = 0;
		EnableMenuItem( item_out );
		EnableMenuItem( item_select );
		EnableMenuItem( item_sort );
		break;
	    case item_out:
		ProcessDROutput( &DR_list, skip_count, DR_type, false, 
		                 theSet.molecule, theSet.struc_count );
		break;
            case item_make:
		if ( OptionTrue( opt_sterlib ) )
		    ModOption( opt_rstinp, "user" );  /* force user libraries; cf. AquaPseudo */
                MakeDRIndices( &DR_list, &theSet );
		if ( GetOptionSubItem( opt_sort ) != sort_off )
		    SortRestraints( &DR_list, GetOptionSubItem( opt_sort ) == sort_class ? true : false );
        
		skip_count = CleanupDRList( &DR_list, skip_count, &exit_flag );
        
		dist_inf = MakeDistinf( &DR_list, skip_count );
		dist_rst = MakeDistrst( &DR_list, skip_count );
		dist_count = DR_list.dist_count;
		ftype = DR_list.struc_count>1 ? inf_multi_model : inf_normal;
		CreateAquaDistFiles( DR_type, ftype, "Opening distance files", true, false );
		WriteDistinfFile( dist_inf, dist_count-skip_count, DR_list.struc_count, DR_ident, theSet.filnam );
		WriteDistrstFile( dist_rst, dist_count-skip_count, DR_ident, theSet.filnam );
		EnableMenuItem( item_stereo );
		EnableMenuItem( item_sort );
                break;
	    case item_sort:
		SortRestraints( &DR_list, GetOptionSubItem( opt_sort ) == sort_class ? true : false );
		break;
	    case item_select:
		skip_count = SelectRestraints( &DR_list, theSet.molecule );
		break;
            case item_stereo:
		oldmenu = SetMenu( &menu_stereo, menulist_stereo );
		err = false;
		while ( !done2 )
		{
		    item2 = ReadMenuItem();
		    switch ( item2 )
		    {
			case item_satur:
			    average = false;
			    satur = true;
			    done2 = true;
			    break;
			case item_assign:
			    average = false;
			    satur = false;
			    done2 = true;
			    break;
			case item_average:
			case item_sum:
			    average = true;
			    satur = false;
			    done2 = true;
			    break;
			case item_QUIT:
			case item_NULL:
			    err = true;
			    done2 = true;
			    break;
			default:
			    ProcessMenuError( menu_fatal );
			    break;
		    }
		}
		ResetMenu( oldmenu );
		if ( err )
		    break;
		dist_infl = MakeDistinfl( average, satur, dist_inf, dist_count-skip_count, &theSet, &a_count );
		if ( !dist_infl )
		    break;
		strcpy( dinf_filnamNormal, dinf_filnamG ); /* save the name of the original (normal) file */
		CreateAquaDistFiles( DR_type, inf_extended, "Opening extended distance files", true, false );
		WriteDistrstFile( dist_rst, dist_count-skip_count, DR_ident, theSet.filnam );
                WriteDistinflFile( dist_infl, dist_count-skip_count, a_count, 
		                   DR_ident, theSet.filnam );
		break;
            case item_show:
                WriteStrucset( stdout, &theSet );
                break;
            case item_QUIT:
            case item_NULL:
                exit( exit_flag );
                break;
            default:
                ProcessMenuError( menu_fatal );
                break;
        }    
    }
    exit( exit_flag );
}

/****************  MakeDRIndices  ****************/
/*
**  Find the atom indices of the atoms in the distance restraint list.
**  The routine also sets the new sequential residue numbers of these atoms.
**
**  If an atom is not found the atom index is set to -1, and if <opt_skip>
**  has not been set the program exits.
**
**  Name guessing - as implemented in FindAtom - is tried if <opt_guess>
**  has been set.
*/

void MakeDRIndices( DRList *DR_list, Strucset *setptr )
{
    char       *fct_name = "MakeDRIndices";
    DRestraint *DRs;
    int        struc_num, DR_count;
    Mol        *molptr = setptr->molecule;
    Atom       *at;
    Boolean    guess, skip, err = false, stop = false;
    int        find_which, print_which;


    skip = OptionTrue( opt_skip );
    guess = OptionTrue( opt_guess );
    
    if ( terseness < 6 ) {
        STATUS( "finding the restrained atoms in the atom list" );
    }
    if ( OptionTrue( opt_oldmap ) )
    {
	find_which = use_ext;
	print_which = use_ext;
        if ( terseness < 6 ) {
	    STATUS( "using the external atom names" );
        }
    }
    else
    {
	find_which = use_int;
	print_which = use_both;
        if ( terseness < 6 ) {
	    STATUS( "using the internal atom names" );
        }
    }
    if ( DR_list->struc_count > 1 )
        if ( terseness < 6 ) {
	    STATUS( "detailed report only for first restraint model" );
        }

/* loop over the models */

    DRs = DR_list->list;
    for ( struc_num = 1; struc_num <= DR_list->struc_count; struc_num++ )
    {

/* loop over the restraints */

	DR_count = DR_list->dist_count;
	while ( DR_count-- )
	{
	    if ( (at = FindAtom2( &DRs->idx_A, &DRs->atom_A, molptr, find_which, guess )) )
		DRs->atom_A.res_num_new = at->res_num_new;
	    else
	    {
		if ( struc_num == 1 )
		    WriteNotFoundAtom( stdout, &DRs->atom_A, find_which, print_which );
		DRs->idx_A = -1;
	    }
    
	    if ( (at = FindAtom2( &DRs->idx_B, &DRs->atom_B, molptr, find_which, guess )) )
		DRs->atom_B.res_num_new = at->res_num_new;
	    else
	    {
		if ( struc_num == 1 )
		    WriteNotFoundAtom( stdout, &DRs->atom_B, find_which, print_which );
		DRs->idx_B = -1;
	    }
    
	    if ( DRs->idx_A == -1 || DRs->idx_B == -1 )
	    {
		if ( skip )
		    err = true;
		else
		{
		    stop = true;
		    break;
		}
	    }
	    DRs++;
	}
	if ( err || stop )
	{
	    printf( " * One or more restraint atoms could not be found" );
	    if ( DR_list->struc_count > 1 )
		printf( " (restraint model %i)\n", struc_num );
	    else
		putchar( '\n' );
	}
	if ( stop )
	    exit( 1 );
    }
}

/****************  MakeDistinf  ****************/
/*
**  Create a Distinf list from a DRList, which has been
**  processed by MakeDRIndices.
**  All restraint models are processed, concatenating the results
**  into a long list (cf. format of Distval lists).
*/

Distinf *MakeDistinf( DRList *DR_list, const int skip_count )
{
    char        *fct_name = "MakeDistinf";
    int         count;
    DRestraint  *DRs;
    Distinf     *dist_inf, *d_inf;

    dist_inf = CREATE_NEW( Distinf, DR_list->struc_count * ( DR_list->dist_count - skip_count ) );
    if ( !dist_inf )
	NOMEM( "dist_inf" );
    
    d_inf = dist_inf;
    DRs = DR_list->list;
    count = DR_list->struc_count * DR_list->dist_count;
    while ( count-- )
    {
	if ( DRs->keep )
	{
            d_inf->dtype = 0;
	    d_inf->idx_A = DRs->idx_A;
	    d_inf->idx_B = DRs->idx_B;
	    d_inf->res_A = DRs->atom_A.res_num_new;
	    d_inf->res_B = DRs->atom_B.res_num_new;
	    d_inf++;
	}
        DRs++;
    }

    return( dist_inf );
}

/****************  MakeDistrst  ****************/
/*
**  Create a Distrst list from a DRList, which has been
**  processed by MakeDRIndices.
**  Only the first restraint model is processed, since all
**  models should have the same list. 
*/

Distrst *MakeDistrst( DRList *DR_list, const int skip_count )
{
    char        *fct_name = "MakeDistrst";
    int         count;
    DRestraint  *DRs;
    Distrst     *dist_rst, *d_rst;

    dist_rst = CREATE_NEW( Distrst, DR_list->dist_count - skip_count );
    if ( !dist_rst )
	NOMEM( "dist_rst" );
    
    d_rst = dist_rst;
    DRs = DR_list->list;
    count = DR_list->dist_count;
    while ( count-- )
    {
	if ( DRs->keep )
	{
	    switch ( DRs->stat )
	    {
		case bnd_upper:
		    d_rst->typ_bnd = 'u';
		    break;
		case bnd_lower:
		    d_rst->typ_bnd = 'l';
		    break;
		case bnd_uplo:
		    d_rst->typ_bnd = 'b';
		    break;
		case bnd_uld:
		    d_rst->typ_bnd = 'B';
		    break;
		default:
		    FATAL( "Program error - illegal value DRs->stat" );
	    }
            d_rst->lo_bnd = DRs->lo_bnd;
	    d_rst->up_bnd = DRs->up_bnd;
	    d_rst->dist   = DRs->dist;
	    d_rst++;
	}
        DRs++;
    }

    return( dist_rst );
}

/****************  MakeDistinfl  ****************/
/*
**  Create an extended Distinf list from a normal Distinf list.
*/

#define max_FORM   6

Distinfl *MakeDistinfl( const Boolean average, const Boolean satur, Distinf *old_list, 
                        const int dist_count, Strucset *set, int *a_count )
{
    char      *fct_name = "MakeDistinfl";
    char      *lib_nam_temp;
    Distinfl  *dist_inf, *dinf;
    Distidx   *dist_idx, *didx, *dx;
    Mol       *molptr;
    Atom      *atoms, at;
    Psatom    *stereo_list, *stereo;
    int       atom_count, stereo_count, err_count = 0, count, i, max, idx, *eq;
    int       **equiv_ptrlist, **equiv_ptr, *type_list;
    int       *equiv_list, *equiv;

    molptr = set->molecule;
    atom_count = molptr->atom_count;

/* generate the required "stereo" information about equivalent atoms */

    if ( average ) {
        lib_nam_temp = "PseudoLIB-aqua";
    } else {
        lib_nam_temp = "StereoLIB-aqua";
    }
    stereo_list = ProcessStereoLIB( &stereo_count, molptr, lib_nam_temp);
    if ( !stereo_list || stereo_count == 0 )
    {
	puts( " * No stereo definitions!" );
	return( NULL );
    }

/* create a list of indices to equivalent atoms */

    max = max_FORM + 2;         /* one extra for pseudo, one for terminator */
    equiv_list = CREATE_NEW( int, max * stereo_count );
    if ( !equiv_list )
	NOMEM( "equiv_list" );

/* create a list of pointers to elements of the previous list */
/* and a list of types */
/* both lists define atom properties */

    equiv_ptrlist = CREATE_NEW( int *, atom_count );
    if ( !equiv_ptrlist )
	NOMEM( "equiv_ptrlist" );
    type_list = CREATE_NEW( int, atom_count );
    if ( !type_list )
	NOMEM( "type_list" );

/* loop over the stereo definitions (cf. MakePseudoSet) */
/* generate and store the atom indices */

    equiv = equiv_list;
    stereo = stereo_list;
    count = stereo_count;
    while ( count-- )
    {
	eq = equiv;
	if ( stereo->num_form > max_FORM )
	{
	    sprintf( msg_string, "too many constituent atoms (%i)", stereo->num_form );
	    ERROR( "change max_FORM" );
	    return( NULL );
	}
	for ( i = 0; i < stereo->num_form; i++ )
	{
/* check whether constituent atom can be found */
	    at = FillAtom( stereo->res_id, stereo->res_nam, 
			   stereo->form_nam[i], use_int );

	    if ( !FindAtom2( &idx, &at, molptr, use_int, false ) )
	    {
		WriteNotFoundAtom( stdout, &at, use_int, use_int );
    
		at = FillAtom( stereo->res_id, stereo->res_nam,
			       stereo->pseudo_nam, use_int );
		printf( "Cannot apply stereo definition for " );
		WriteAtom( stdout, &at, use_int );
		putchar( '\n' );
		err_count++;
		*equiv = -1;       /* disable this definition */
		goto next;
	    }
/* store the index of this equivalent atom */
	    *eq++ = idx;
	}
/* now add the pseudo atom itself */
	at = FillAtom( stereo->res_id, stereo->res_nam,
		       stereo->pseudo_nam, use_int );
	if ( !FindAtom2( &idx, &at, molptr, use_int, false ) )
	{
	    WriteNotFoundAtom( stdout, &at, use_int, use_int );
	    puts( "Cannot apply stereo definition for this atom" );
	    err_count++;
	    *equiv = -1;           /* disable this definition */
	    goto next;
	}
	*eq++ = idx;
/* terminate the definition */
	*eq = -1;

/* store pointer to this definition and store the pseudo types */
/* for all atoms involved except the pseudo atom itself, or only for the pseudo */
/* (depending on later usage: for Assign or for DoAverage) */
	eq = equiv;
	while ( *eq != -1 )        /* just to be sure to stop in time */
	{
	                           /* (don't) store (only) for the last i.e. the pseudo itself */
	    if ( ( *(eq+1) != -1 && !average ) || ( *(eq+1) == -1 && average ) || satur )
	    {
		equiv_ptrlist[ *eq ] = equiv;
		type_list[ *eq ] = stereo->pseudo_typ;
	    }
	    eq++;
	}
/***
	{
	    Boolean done=false;
	    eq = equiv;
	    while ( !done )
	    {
		done = (*eq == -1);
		printf( " %i", *eq++ );
	    }
	    putchar( '\n' );
	}
***/
next:
	equiv+=max;
        stereo++;
    }

/* report results */

    if ( err_count )
    {
        if ( err_count == stereo_count )
	{
    	    printf( " * No stereo definitions could be made\n" );
	    exit_flag = 30;
	    return( NULL );
	}
        else
	{
	    printf( " * %i stereo definitions could not be made\n", err_count );
	    printf( " * %i stereo definitions have been made\n", stereo_count - err_count );
	    if ( !average ) 
		exit_flag = 40;
/* 
** NOTE: this is only an error when the undefined pseudo occurs in the 
** restraint list; the latter will not be the case, since the normal distinf 
** list cannot contain references to unexpandable pseudo atoms.
** When average==false the situation is different, since the pseudo atom
** need not be present in the original restraint list at all.
*/
	}
    }

/* print the list (checking only) */
    if (terseness < 6 ) {
        atoms = molptr->atoms;
        equiv_ptr = equiv_ptrlist;
        count = atom_count;
        while ( count-- )
        {
	    if ( *equiv_ptr )
	    {
	        WriteFAtom( stdout, atoms, use_both );
	        eq = *equiv_ptr;
	        while ( *eq != -1 )
	        {
		    if ( average && *(eq+1) == -1 )  /* skip the pseudo? */
		    {
		        eq++;
		        continue;
		    }
		    printf( " %i", *eq++ );
	        }
	        putchar( ' ' );
	        eq = *equiv_ptr;
	        while ( *eq != -1 )
	        {
		    if ( average && *(eq+1) == -1 )  /* skip the pseudo? */
		    {
		        eq++;
		        continue;
		    }
		    PRINTIDv( molptr->atoms[ *eq ].atom_int );
		    eq++;
	        }
	        putchar( '\n' );
	    }
	    equiv_ptr++;
	    atoms++;
        }
    }

/* allocate space for the Distinfl list */

    dist_inf = CREATE_NEW( Distinfl, dist_count );
    if ( !dist_inf )
	NOMEM( "dist_inf" );

/* allocate space for the Distidx data */

    dist_idx = CREATE_NEW( Distidx, dist_count * 2 * (max - 1) );
                                                /* don't need to store terminator */
    if ( !dist_idx )
	NOMEM( "dist_idx" );

/* now make the Distinfl list */

    *a_count = 0;
    dinf = dist_inf;
    didx = dist_idx;
    count = dist_count;
    while ( count-- )
    {
	dinf->dtype = 0;

	dx = didx;
	dinf->A = dx;

	idx = old_list->idx_A;
	eq = equiv_ptrlist[ idx ];
	if ( eq )
	{
	    i = 0;
	    dinf->ass_A = -1;
	    while ( *eq != -1 )
	    {
		if ( average && *(eq+1) == -1 )  /* skip the pseudo? */
		{
		    eq++;
		    continue;
		}
		if ( *eq == idx )
		    dinf->ass_A = i;
		dx->res = old_list->res_A;
		(dx++)->idx = *eq++;
		i++;
	    }
	    dinf->count_A = i;
	    dinf->atype_A = type_list[ idx ];
	}
	else
	{
	    dinf->ass_A = 0;
	    dinf->count_A = 1;
	    dinf->atype_A = 0;
	    dx->res = old_list->res_A;
	    (dx++)->idx = idx;
	}

	dinf->B = dx;

	idx = old_list->idx_B;
	eq = equiv_ptrlist[ idx ];
	if ( eq )
	{
	    i = 0;
	    dinf->ass_B = -1;
	    while ( *eq != -1 )
	    {
		if ( average && *(eq+1) == -1 )  /* skip the pseudo? */
		{
		    eq++;
		    continue;
		}
		if ( *eq == idx )
		    dinf->ass_B = i;
		dx->res = old_list->res_B;
		(dx++)->idx = *eq++;
		i++;
	    }
	    dinf->count_B = i;
	    dinf->atype_B = type_list[ idx ];
	}
	else
	{
	    dinf->ass_B = 0;
	    dinf->count_B = 1;
	    dinf->atype_B = 0;
	    dx->res = old_list->res_B;
	    (dx++)->idx = idx;
	}

	*a_count += dinf->count_A + dinf->count_B;
	didx += 2 * (max - 1);
	old_list++;
	dinf++;
    }

/* clean up */

    free( type_list );
    free( equiv_ptrlist );
    free( equiv_list );

    return( dist_inf );
}

